Weakly polydisperse systems: Perturbative phase diagrams that include the critical 

region 
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The phase behaviour of a weakly polydisperse system, such as a colloid with a small spread of 
particle sizes, can be related perturbatively to that of its monodisperse counterpart. I show how this 
approach can be generalized to remain well-behaved near critical points, avoiding the divergences 
of existing methods and giving access to some of the key qualitative features of polydisperse phase 
equilibria. The analysis explains also why in purely size-polydisperse systems the critical point is, 
unusually, located very near the maximum of the cloud and shadow curves. 
PACS numbers: 64.10.+h, 68.35.Rh, 82.70.-y, 05.20.-y 



Most soft materials are inherently polydisperse: spher- 
ical colloids have a spread of diameters, polymers an 
effectively continuous distribution of chain lengths etc. 
The effects of polydispersity on equilibrium phase be- 
haviour are both of fundamental interest and of rele- 
vance to production, processing and shelf life of food- 
stuffs, personal hygiene products, paints and coatings etc. 
Careful sample preparation can reduce polydispersity, to 
e.g. diameter standard deviations of around 5% for col- 
loids [X], but not eliminate it. Existing theories that cal- 
culate the effects of such weak polydispersity by expand- 
ing around the phase equilibria of a monodisperse ref- 
erence system 2, 3] have proved useful, but break down 
around critical points (CPs) where the relevant prefactors 
scale with the compressibility. The critical region can be 
studied alternatively by a Landau expansion [4]; while 
not limited to weak polydispersity, this is impractical for 
all but the simplest models. Generic tools for predict- 
ing polydispersity effects also in the critical region would 
clearly be useful: as we will see, the qualitatively dis- 
tinct features of polydisperse phase behaviour manifest 
themselves most clearly here, and critical parameters can 
be estimated reliably from simulations 0, Q and used to 
constrain theoretical models. My aim here is to supply 
such a tool, by formulating a perturbative approach that 
remains applicable in the critical region. The approach 
allows for the prediction of the full phase behaviour of 
generic weakly polydisperse materials from that of the 
monodisperse reference system. As a by-product, the 
analysis rationalizes previously unexplained peculiarities 
in the phase behaviour of purely size-polydisperse sys- 
tem. The key idea is that, because a polydisperse system 
has a critical temperature different from its monodisperse 
counterpart, one has to allow perturbations not just of 
the coexisting densities but also of temperature. 

Denote the polydisperse attribute, e.g. particle diam- 
eter, by (J, and its average in the parent phase that is 
allowed to phase separate by ctq. In terms of the scaled 
deviations from the parental mean, e = (cr — cro)/<''0: the 
normalized parent size distribution p{e) then has zero 
mean and, by the assumption of weak polydispersity, 
small variance w — Jdep{e)e^. We describe a generic 



phase of the system by its (number) density distribution 
p(e). Its excess Helmholtz free energy density in units 
of fceT, /poly, is a functional of p{e) and an ordinary 
function of inverse temperature (3 = 1/k^T. For weak 
polydispersity this can be expanded up to 0{e^) as 0,0]: 



/poiy([p(e)], /3) = / + Ap, + Bp2 + Cp\ 



(1) 



Here the monodisperse reference free energy / is a func- 
tion of the overall number density pq — Jde p{e) and /3, 
as are the coefficients A, B, C; we abbreviate p = po 
below. The higher-order moment densities are defined 
as pn — Jdep{£)e^ {n — 1,2). From /poly one finds the 
excess chemical potentials jl{t) = Sfpo\y/Sp{e) 

= fp + ApPi + BpP2 + Cppi + {A + 2Cpi)e + Be^ 

where subscripts indicate derivatives, i.e. fp = df /dp 
etc. Equality of the full chemical potentials p{e) — 
Inp(e) -|- /2(e) then implies that phases in coexistence 
have density distributions differing by a factor of the 
form cxp(Ao + Aie + A2e^), as expected from the mo- 
ment structure of ([1]) 0, [§]■ The pressure Hpoiy — 
p+ Jde p{e)p{e) - /poiy becomes 



Hpoiy = n + Apppi + Bppp2 + (Cpp + C)pI 



(2) 



with the monodisperse reference pressure IT = p+pfp — f- 
In a polydisperse system the monodisperse binodal 
splits into cloud and shadow curves The former 

records the density of the majority phase at the onset of 
coexistence, the latter that of the incipient phase. The 
cloud phase is then identical to the parent, with a density 
distribution of the form pp(e), while that of the coexist- 
ing shadow can be written pp(e)exp(Ao + Aie -I- A2e^). 
Equahty of pressures and of the 0(e"), O(e^) and 0{e^) 
terms in the chemical potentials p(e) gives the conditions 



n + wBpp"^ 
fp + wBpp 
A 
B 



n' + A\^p'p[ + B'^p'p'^ + (C>' + C")p'i^ 

^'p + A'^p'^+B'^p'2 + C'pP'^' + \^ 
A' + 2C'p'i + Ai 
B' + \2 
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where unprimed/primed quantities relate to the cloud 
and shadow phase, respectively, e.g. p', p^, p'^ are den- 
sity moments in the shadow and Ap = Ap(/9, /?), = 
Ap{p',l3); we have used that pi = p{e) = and p2 = 
p{e^) — wp in the cloud. 

The basis of the perturbation expansion is the 
monodisperse phase coexistence, with coexisting densi- 
ties p, p' = pexp(Ao) at temperature /3. One then 
has w ~ p'l = p'2 = Q and the Lagrange multipliers 
become \o ^ fp - f'p = Afp, Xi = A - A' = AA, 
A2 = B — B' = AB. The A-notation indicates differences 
between coexisting monodisperse phases. When the par- 
ent is made polydisperse, p is perturbed to the cloud 
point density pc = p + Spc, and the Lagrange multipliers 
to Ao = Afp + SXo, Ai = AA + ^Ai, A2 = AB + S\2. Ex- 
panding to linear order in all small quantities including 
w, the shadow moments become p'^ — p' + Sp'^ with 

5p'^ = p'{5pjp + 6X0 + w[^{AAf + AB]} (3) 

and p'l = wp'AA, p'2 — wp' . One now inserts these into 
the four phase coexistence conditions above and expands, 
crucially also allowing a temperature shift /3 (3 + 6f3. 
The fourth conditions only serves to determine 6X2 so 
can be omitted. The third condition reduces to 

ApSpc - A'pdp', + AAf) 6(3 = 2wp'C'AA + SXi (4) 

while the first two yield, after eliminating SXq in favour 
of Sp's using ([3]), 

np6p, - Il'pSp'^ + Anp 5(3 = 

= w[App'^AA - AiBpp")] 
ifpp + l/p)Sp, - Cf'pp + ^/P')6p', + Afpp 6(3 = 

= w{A'pp'AA - {AAf/2 - A{Bpp) - AB} 

The cofficients involving / and H can be written in terms 
of the full monodisperse free energy density including the 
ideal part, / = + Inp) + f: fpp + l/p = fpp, fpp = 
fp/3, Up = pfpp and 11/3 = pfpf3 - //3- Taking appropriate 
linear combinations to eliminate 5pc or 6 pi. then gives 

fppSpc + [fpf3 - Afp/Ap] 5(3 - (5) 

= w{p'[\{AAf + AB]IAp-pBp} 
f'pp5p', + [f'pp^Afp/Ap]5(3^ (6) 

= w{p[\{AAf + AB]/Ap - p'B'p - p'A'pAA} 

If temperature is held constant ((5/3 = 0), one retrieves 
Evans' results [3] for the cloud and shadow density shifts. 
This approach breaks down at the CP because 1 / fpp and 
1/ fpp, which are proportional to the compressibilities in 
the monodisperse reference phases, diverge. 

In the more general setup here we can choose a pre- 
scription for SXi — 0{w) and deduce the temperature 
shift 5(3, which will be of the same order. This freedom 
does not mean that the perturbation theory is ill-defined. 



In fact, at a fixed perturbed temperature (3^ = (3 5(3 
away from the CP, the results are independent of 5(3. To 
see this for e.g. the cloud point density, write ([5]) as 

Pe(/3o) = p(/3)+P/3(/3)5/3 + u;i?(/3) (7) 

The function wR{(3) is the r.h.s. of (O divided by fpp, 
with the /3-dependence of the monodisperse coexisting 
densities p and p' inserted. That the prefactor of 5(3 must 
equal the derivative pp = dp/d(3 along the monodisperse 
binodal follows by setting w = 0, in which case there is 
no additional polydispersity-induced shift. Now if 5(3 = 
[3], pc(/3o) = p{Po) + wR{(3o). This agrees to 0{w) 
with (O, for any temperature shift 5(3 = 0{w), since 
p{Po) = p{(3) + Pi3{(3)5(3 + 0{w'^); the difference between 
R{(3) and R{Pq) likewise only gives O(w^) corrections. 

To guide the choice of 5Xi one can use the intu- 
ition that in a reliable perturbation theory correspond- 
ing monodisperse and polydisperse state points should 
be mapped to each other. At the CP, where cloud and 
shadow become identical, also their size distributions do, 
so we must have Ai = AA -I- <5Ai = i5Ai = 0. This sug- 
gests taking 5Xi = throughout the phase diagram; to 
eliminate the function C from ([4]), however, we choose 

5Xi = -2wp'C'AA (8) 

which still vanishes at the CP. 

Equations with the choice ([5]) constitute our gen- 
eralized perturbation theory; the set of three linear equa- 
tions is easily solved for 5pc, 5p'^, 5(3. I show below that 
the resulting perturbations do remain finite in the crit- 
ical region for systems with classical critical behaviour. 
This encompasses a large variety of models used to de- 
scribe polydisperse complex fluids, such as van der Waals 
theory for liquids, free volume theory for colloid-polymer 
mixtures, the Flory-Huggins and Wertheim's statistical 
associating fluid theories for polymers etc (see references 
in [9]). For monodisperse reference systems with non- 
classical critical behaviour one expects that polydisper- 
sity will Fisher- renormalize the critical exponents [lO|, 
multiplying them by 1/(1 — a). For small w a crossover 
from the monodisperse to the polydisperse exponents 
should then occur very close to the CP, and this will 
not be accounted for by our perturbation theory. Study- 
ing this crossover is an interesting issue for future work. 
Because the specific heat exponent a « 0.11 is small in 
d = 3, however, the quantitative effects of Fisher renor- 
malization may well be small enough for the perturbative 
approach to remain an accurate approximation. 

To obtain an explicit expression for the CP shift due to 
polydispersity, we start by assuming a smooth (Landau- 
like) expansion of / near criticality. The coexisting den- 
sities then go as p - = -{p' - p*) = [2{(3 - (3*)/bY/^ 
to leading order, with b some constant; asterisks denote 
quantities at the CP. Defining r = [2{(3 - (3*)/bY/'^, we 
can then write p — p* = r -\- ar^, p' — p* — — r -|- or^. 



3 



/? = /?* + ^br^ to ©(r^). The coefficient b is given by 
b = —fpppp/ (3/*p^); a will cancel to the order we require. 
To expand our equations for small r write ([4]), together 
with (|5|6p divided by Ap, in matrix form MS — wv. 
Here = {6pc, Sp'^, S(3) and d gathers the appropriate 
right-hand sides. Substituting the above expansions for 
p, p', (i and using f*p = f*pp = at the CP then yields 
M = M* + M;r + . . ., V = V* + v*r + . . . with 



M* 





(9) 



where V = \p*{A*p^ - B*^) - B*. We want to deduce 
from this the expansion of the solution vector, S ^ 5* + 
rS* + . . . Because M* is degenerate, the zeroth order 
M*5* — wv* fixes the 5/3-component of 5* as 



6p* = w[p*{Af-B;p)-2B;]if;pp 



(10) 



and otherwise (only) tells us that the other components 
must be equal, 5pc — Sp'^ = 6p* . This is as it should be: 
we built the perturbation theory to map the monodis- 
perse CP to the polydisperse CP, where pc = p'^- Math- 
ematically it follows from the vanishing of the first com- 
ponent of V* , confirming the intuition that SXi needs to 
be chosen to make the r.h.s. of ((4]) zero at the CP. 

To determine Sp*, we use the 0{r) condition M*S* + 
M;S* = wv;. Multiplying from the left by = (0, 1, 1) 
eliminates the first term and gives 6"""^^*^* = we^v^. 
By inserting the explicit expressions for M* and from 
the small-r expansion and the form 6* — {Sp*, 6p*,6(3*)'^ 
one then finds 

f;ppp5p*+rpppp5p* = w[iAf~zB;p+p*{iA;A;p-B;j] 

Together with PH)) this gives the desired shifts of the 
critical density, 5p* , and critical temperature, 513* \ these 
results can be checked to agree with the general poly- 
disperse CP criterion Q, expanded directly for small w. 
Both shifts 5p* and Sfi* are ©(w), so that the loci of the 
CPs for increasing w depart from the monodisperse CP 
along a line of nontrivial slope in the [p, (3) plane. 

Taking the expansion further, we can obtain the slope 
of the cloud and shadow curves at the CP, and the lo- 
cation of their maxima. Both are related to (5/3* — 
d((5/3)/dr|r=o, the (5/3-component of This follows af- 
ter some algebra from the first order equation M*S; + 
M*S* = wv; together with the known S* as 



Sf3: = -wA;i2A; + p*A;p)/f;pp 



(11) 



The slope of the cloud curve, plotted as /3 vs pc, is 
d{P + SP)/d{p + Spc) = [diP + SP)/dr][d{p + 6pc)/dr]-\ 
At the critical point, the first factor equals i5/3*. This is 
C(w), so we only need the Oiw'^) contribution, dp/dr — 
1, of the second factor to find the cloud slope at the CP 



as 6p; . The shadow slope has the same modulus but op- 
posite sign since dp' /dr = — 1. For the cloud and shadow 
maxima we want to find r such that d{f3 + 6(3) /dr — 0, 
giving br + 6f3; = oi r — —5(3;/b. Because this is 0{w), 
it is correct to replace 5(3r by (5/3* as we have done. It 
follows that the cloud and shadow maxima are to 0{w) 
located at pc,s = p* + Sp* =F Sf3;/b, i.e. equidistant either 
side of the polydisperse CP. 

Finally we turn to purely size-polydisperse ("scal- 
able" [s*]) systems, in which the interaction potential re- 
mains unchanged when the diameters ai and position 
vectors of all particles are scaled by a common factor. 
For pair interactions, this holds if the dependence on dis- 
tance Tij of the potential between two particles i and j 



can be written as Vij (r^. 



^{nj/£{ai,a.j)) with the 



interaction distance £ a homogeneous function of first 



degree like £= {ai + aj)/2 or i 



This 



eludes e.g. polydisperse (additive) hard spheres, or size- 
polydisperse Lennard- Jones mixtures. Simulations and 
theoretical calculations 0, |^ [ill, [l3] show that for such 
systems, unusually in the presence of polydispersity, the 
critical point is located very near the maximum of the 
cloud and shadow curves (and thus also close to the maxi- 
mum of a plot of pressure versus temperature at the onset 
of phase coexistence [i,[iH). Moreover, plotted in terms 
of volume fraction instead of density, cloud and shadow 
curves almost coincide d, H, [ll, 12 1. 

This behaviour is fully predicted by our 0{w) pertur- 
bation theory. As regards the location of the CP, we need 
only use that in size-polydisperse systems [3] 



A = 3{U/p-l) = 3{fp-f/p-l) 



(12) 



Working out Ap and App, substituting into pT|) and using 
f *p = fppp then gives (5/3* = 0. The cloud and 
shadow slopes at the CP, ±(5/3*, thus vanish to 0{w), so 
that the CP must be at their maximum as claimed. 

To see the coincidence of the volume fraction repre- 
sentations of cloud and shadow, we measure densities in 
units of (7r(7Q/6)^^ so that the volume fraction can be 
written (p = Jdep{e){l + e)'^. Expanding to 0{w) gives 
for the cloud volume fraction (jjc = p(l + Sw) -I- Spc, while 
for the shadow (j)s — p'{l + 3wAA + 3w) + 6p'^. One can 
now write equations (jlHS]) in terms of the volume frac- 
tion shifts (5(/>c = 4>c — p and (5(/)g = (f)'^ ~ p' . Swapping 
p and p' gives equations for 5<j)'^ and 5(j)s', coincidence 
requires that 5(j)c — S(f>s and 6(j)'^ — S(f)i,, i.e. any point 
on the monodisperse binodal is perturbed into cloud and 
shadow points with identical volume fractions. The orig- 
inal and swapped equations turn out to differ by a term 
proportional to pAA/Ap + ip'f'pp — p'A'p on the r.h.s.; 
from (fT2)l this indeed vanishes for size-polydisperse sys- 
tems. To ensure that also the temperature shift (5/3 is the 
same, one needs to modify slightly the choice ^ of (5Ai, 
taking e.g. 



(5Ai = wp'AA{~2C' + 3A') 



(13) 
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FIG. 1: (a) Cloud and shadow curves for a model of a 
Lennard- Jones mixture with polydispersity in both parti- 
cle diameters and interaction strengths. Symbols: full non- 
perturbative solution. Thick lines: perturbation theory with 
choice (|13[) : using (|8]) only gives noticeable differences for 
the high-density part of the shadow curve (dash-dotted line) . 
Thin lines: constant-/? approach of [3] • Dotted line: monodis- 
perse binodal. (b) Similar plot for a purely size-polydisperse 
mixture, in the volume fraction representation. Dash-dotted 
line: unregularized perturbation theory. 



sentation cloud and shadow curves coincide as predicted, 
and the CP is very near their maximum. Again, our 
approach avoids divergences in this region and is quanti- 
tatively accurate even for the substantial degree of poly- 
dispersity w^/^ — 0.2 considered here, i.e. a size standard 
that is 20% of the mean. The dash-dotted line gives the 
raw results of the perturbation theory; the unphysical 
upward curvature comes from an increase in the temper- 
ature shifts \5I3\. This illustrates that, away from the CP 
where we should have 5\i = 0, it is not obvious how to 
choose 5\i to extract the most reliable predictions for 
finite w. Intuitively, one expects temperature shifts to 
be needed primarily in the critical region and to become 
smaller elsewhere. This can be enforced by e.g. taking the 
raw 5(i and using a soft thresholding function so that val- 
ues with I (5/3 1 > \5(3*\ are reassigned to be no larger than 
2 1 5/3* I; 5pc and 5p'^ are then recalculated from (|5I6[) . The 
thick line in Fig. [1] shows that this gives accurate and 
physically sensible results. 

The approach presented here allows the calculation of 
full phase diagrams for generic weakly polydisperse sys- 
tems, giving access to key features such as the separa- 
tion of the maxima of cloud and shadow from the CP, 
the differing cloud and shadow slopes at the CP, and the 
polydispersity-induced shifts of the critical parameters. 
Even for excess free energies with moment structure it 
can be helpful in avoiding the computational complexi- 
ties of a full phase equilibrium calculation while re- 
maining quantitatively accurate in a substantial range of 
polydispcrsities (up to 20% in Fig. [J) that includes typ- 
ical experimental values [ij. The theory also generalizes 
easily to phase equilibria inside the coexistence region, 
as well as systems with several polydisperse attributes. 



To illustrate the theory we show in Fig.[TJa) the predic- 
tions for a Lennard-Jones mixture with polydispersity in 
both particle diameters and interaction strengths, mod- 
elled with a moment excess free energy (tIIsI as in [l3. 14|. 
for a uniform ("top hat") size distribution p(e) of width 
■u;i/2 = 0.07. The generalized perturbation theory pro- 
duces cloud and shadow curves with smoothly connecting 
sub- and super-critical branches, avoiding the divergences 
near the monodisperse CP of the constant-/? approach 
Comparing with the symbols, our perturbation theory is 
quantitatively closer to the full solution of the polydis- 
perse phase equilibrium conditions also elsewhere, in par- 
ticular for the high-density branch of the shadow curve. 
The overall level of agreement is quite remarkable given 
that the shifts from the monodisperse reference system 
(dotted line) are substantial. 

Fig. [TJb) shows similar results for a mixture with pure 
size polydispersity @, [3| ■ In the volume fraction repre- 
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